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We consider a system composed of two identical light quarks (qq) and two identical antiquarks 
(QQ) that can be linked either as two mesons or as a tetraquark, incorporating quantum correlations 
between identical particles and an effective many-body potential between particles. We perform a 
3-D Monte Carlo simulation of the system, considering the conhgurations allowed to form: i) Only 
two mesons, ii) Only tetraquark and hi) two mesons and tetraquark . We characterize each case 
and determine whether it is energetically more favorable to form a tetraquark or two mesons, 
as a function of the interparticle separation distance which, for a fixed number of particles, can 
be identihed as a particle density. We determine how the two mesons, which dominate the low 
density regime, mixes with a tetraquark state as the density increases. Properties like the mean 
square radius and the two-particle correlation function are found to reflect such transition, and we 
provide a parameterization of the diquark correlation function in the isolated case. We track the 
dynamical flipping among configurations to determine the recombination probability, exhibiting the 
importance of the tetraquark state. We analize the four-body potential evolution and show that its 
linear behavior is preserved, although the slope can reflect the presence of a mixed state. Results 
are shown for several light-quarks to heavy-antiquarks mass ratios whenever they are found to be 
relevant. 

PACS numbers: 14.40.Rt, 12.38.Lg, 12.39.Pn, 12.39.Jh 


I. INTRODUCTION 

The formation of multi-quark systems can have impor¬ 
tant implications in the phenomena we observe in nature, 
from an enhanced spectroscopy to quark recombination 
effects. A tetraquark state, the simplest of this kind of 
systems, has been extensively studied from the theoreti¬ 
cal and experimental point of view. Recently, experimen¬ 
tal research has provided strong evidence on the forma¬ 
tion of such state [in. These findings signal a new era in 
hadron spectroscopy, where more refined modeling and a 
better understanding of the strong interaction, in the low 
energy regime, is required. Since the early years of the 
quark model, theoretical studies have been performed to 
inquiry on the existence and stability of the tetraquark 
as an isolated object and how its mixing with a 

meson state can help us to understand the observed spec¬ 
troscopy of states like the a meson pm. Less attention 
has been paid to the features of the tetraquark formation 
as two mesons are forced to approach each other as it 
could happen in a meson-meson collision or, on the op¬ 
posite, when the four quarks are produced very close in 
space as in the WW decay, which eventually freeze out 
to two mesons Cl El. At which stage they turn into 
a tetraquark or mixed state? and how their properties 
reflect such modification? are certainly important ques¬ 
tions. 

In the present work we address these questions using an 
effective model (String-Flip Model) to mimic the strong 
interaction among quarks. For definiteness we consider 
two identical light quarks qq (denoted by u) and two iden¬ 
tical antiquarks QQ (denoted by d, s, c and b, depending 
on the mass ratio respect to the light one). The sem¬ 
inal work by Lenz et al. Q laid down a procedure to 
describe such a system in this context and found its gen¬ 


eral properties, as isolated objects, using a harmonic po¬ 
tential interaction. Here, we perform a 3-D Monte Carlo 
simulation of the system including quantum correlations 
between particles to determine whether it is energetically 
more favorable to form a 4-body state (tetraquark) or 
two mesons, as a function of the interparticle separation 
which, for a fixed number of particles, can be identified 
as a particle density. We are interested in s-wave states 
where tensor and spin interaction effects are expected to 
be negligible in the gross features of the properties. In 
section II, we set the basis to identify a quark-antiquark 
state interacting via a linear potential and elaborate on a 
variational approach. In section III, we describe the dy¬ 
namical quark recombination model and generalize the 
variational two-body wavefunction to a four-body sys¬ 
tem. In section IV, we describe the energy evaluation 
from Monte Carlo simulation. In section V, we show the 
results on: A) The determination of the optimal vari¬ 
ational parameter and system energy, B) The hadron 
radial distribution characterization, C) The two-particle 
correlation functions, among the two mesons and among 
quarks. In particular we provide a parameterization of 
the diquark correlation function and the static structure 
factor, D) The dynamical recombination probabilities 
among mesons with/without tetraquark configurations, 
E) We characterize the four body potential, estimating 
its linear behavior and its strength as an effective contact 
interaction. The conclusions are presented in section VI. 


II. QUARK-ANTIQUARK STATE WITH 
LINEAR POTENTIAL 

Let us start with the description of the exact solution 
and the variational approach of a meson state, composed 
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of a quark and an antiquark of mass mi and m 2 respec¬ 
tively. The strong interaction between this pair can be 
represented by a flux tube [B] , which can be effectively de¬ 
scribed by a linear potential V [fq,r 2 ] = k\ri — f 2 \ = kr 
where k is an interaction constant and r = |fi — r 2 | is 
the relative distance between them. Lattice QCD stud¬ 
ies have confirmed the linear behavior of the interaction 
among quarks at large distances [^ . 


A. Exact solution 

For a three-dimensional space, one can analytically 
solve the corresponding Schrodinger equation: 


corresponding to the zeros of the Airy function , which 
then define the eigen-energies 
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and eigen-functions: 

Pnir) = ^Ai . (10) 

The ground state corresponds to the first zero of the Airy 
function, with energy Eq = 2.3381, in mi = m 2 = m = 
k = 1 units. 


'IL 

2 mi 



{ri,r2) 


E (ri,r 2 ) . 


( 1 ) 


B. Meson variational wave function 


Since the potential depends only on the relative distance, 
we can use the center of mass {R) and relative distance 
(r) vectors: 


- miri -I- m2r2 ^ ^ ^ 

R= -;-; r = ri-r2, 

mi + m 2 

and the wave function can be set as: 
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Previous studies have explored how to approach a vari¬ 
ational wave function that can account for the meson en¬ 
ergy [3-flfll. [2lj. Here, we describe, in similar fashion, 
a family of 3-D variational wave functions that approach 
this value and find analytical expressions for the expected 
variational parameter, energy and mean square radius. 
The general variational wave function takes the follow¬ 
ing form: 


FxM = ( 11 ) 


We can neglect the center of mass contribution and focus 
in the relative distance Schrodinger equation. Consider¬ 
ing a spherical symmetry (^= 0 ), the radial solution can 
be set as U{r) = rp{r), such that the radial equation 
becomes: 

where p = mi m 2 /(mi-I-m 2 ) is the reduced mass. 
Switching to the variable 

77 (C) = C/(r) (5) 

The Schroedinger Eqn. (jT]) takes the following form: 


where A is a normalization term, b and A are variational 
parameters and r is the relative coordinate of the two- 
body system. We have found that using b — 1.74716 and 
A = 0.350132, the wave function reproduces the exact 
minimum energy value within a 1% error. Based on this 
analysis, we fix the variational parameter 6 = 3/2 and 
leave the A parameter as the only free parameter. Then, 
the normalized meson variational wave function becomes: 

Analytical results can be obtained for the properties of 
the meson for this form of the variational wave function. 
The value of A minimizing the energy is: 




= 0, 
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whose general solution is a linear combination of the Airy 
functions Ai(C) and Bi(C). The second function is not 
considered as it has a bad behavior at large distances 
(C —>■ 00 , limAi(C) = 0 and limBi(C) = 00 ). Thus, the 
solution is: 


viO = Ai(C) (7) 

The boundary condition at the origin U{0) =0 requires 

(-|( 2 M)'/^) = 0 , ( 8 ) 


corresponding to an energy of 

35/3fcr(|) 

“ 27/3(yt^)l/3> 


(14) 


which is only around 3% higher than the exact solution 
(see Table m. The mean square radius for the meson can 
be obtained following the general definition: 
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m 2 /mi 

Ao 

Eo 

^exact 

{^M)n 

1 

0.4714 

2.3472 

2.33811 

1.50255 

1.44643 

0.5125 

2.2197 

2.21106 

1.38885 

4.6131 

0.6043 

1.9889 

1.98118 

1.52605 

14.0774 

0.6441 

1.9061 

1.8987 

1.73656 


TABLE I: Optimal variational parameter and energy, the ex¬ 
act energy and mean square radius for different mass ratios 
respect to the lightest u mass, mi. 


In the meson case with particles of mass mi and m 2 , it 
becomes 


(^m)o — 


ml + ml 
(mi + m2Y 4 



where the subindex 0 is used everywhere to denote the 
isolated case. We can extend this analysis for quarks of 
different masses, in particular considering the constituent 
masses of the d, s, c and b quarks, we set the ratio of 
their masses respect to the u quark mass, taken as ref¬ 
erence [12 . In Table III we show the corresponding op¬ 
timal variational parameter (Aq), the variational energy 
(i?o), the exact energy (inexact) variational ra¬ 

dius ((r^)o) for the different mass ratios. We denote u 
mass as mi, and m 2 corresponds to effective d, s,c and 
b masses. These results will be used as reference when 
extending to a four-body system. 


III. FOUR-BODY SYSTEM 

Using the previous results, we can now extend our 
study to a four-body system composed of two identical 
M-like quarks (qq) and two identical antiquarks (QQ), 
which can be d,s,c and &-like. To simplify the nota¬ 
tion, we will refer to this four-body system by qQ. At 
very low density, where the distance between quarks is 
large, the model must reproduce a system of two isolated 
mesons. In addition, as the density increases, the inter¬ 
particle separation becomes small allowing the interac¬ 
tion between the two mesons, which can be represented 
by quark exchange or a truly four-body interaction. All 
over, the Pauli blocking among identical particles must 
be enforced, allowing the system to behave like a free 
Fermi system at very h^h densities. Here, we rely on 
the String-Flip Model [ 3 , [ 2 ^, [IJl which has been also 
used to study dense matter systems [H, , to describe 

the dynamical transition among these regimes. 


A. String-flip Potential 

To model such a system, we need to capture the 
properties of the strong interaction. In particular, the 
strong interaction of a many-body system can be well 


approached by considering that the quarks are connected 
by gluon flux tubes Q according to a configuration pro¬ 
ducing the lowest energy state of the whole system. 

Let us denote the position of the particles by fi and r 2 
and the antiparticles by and fq. The way to link the 
four particles of the system consistent with the QCD re¬ 
strictions of color neutrality is two-fold: 
i) The system is composed of two mesons. Then, using 
the linear potential per pairs U(ri,fy) = klfl — fjl = 
kvij, the total potential can be either 

Uni = U(ri,r3)-k U(r2,r4) (16) 

or 

Un 2 = U(ri,r4) -I- U(r 2 ,r 3 ). (17) 

ii) The system is composed of a tetraquark. In this 
case, the shortest path linking the four particles is given 
by the Steiner-tree, which uses two auxiliary vectors (de¬ 
noted by k and 1) to minimize the length of the config¬ 
uration. The vector k links the diquark sub-system and 
the vector I the anti-diquark one. The diquark and anti¬ 
diquark are considered to belong to the 3 and 3 represen¬ 
tations of the color SU(3) symmetry, respectively. Thus, 
the total system form a color singlet state. The interac¬ 
tion among particles in such representations is attractive 
with a strength a half of the quark-antiquark interaction. 
Since there is no direct link between particles but rather 
an interaction via auxiliary vectors, we use such strength 
as the same among all of them. This form of the poten¬ 
tial has been shown to be consistent with the flux tube 
picture Q and found, in the Lattice QCD framework, to 
follow a linear behavior at large distances [l^]. Thus, the 
potential becomes: 

2 4 

U 4 Q = ^ U(fc, r))/2 + ^ V{1, r,)/2 + U(fc, l)/2, (18) 

i=l j=3 

where the 1/2 factor accounts for the decreasing in the 
coupling strength. When all these configurations are al¬ 
lowed, the potential of the system is chosen as the one 
producing the minimum potential energy in a given con¬ 
figuration: 


V = min{Vmi,Vm2, V^q). (19) 


B. Variational Wave function 

Generalizing the meson result to include all the possible 
configurations at different densities, we propose the fol¬ 
lowing variational wave function, which reproduces the 
two isolated mesons at very low densities, allows four- 
body correlations and incorporates the Fermi correlations 
as the density increases: 


'I'a = ^FG (e , 


( 20 ) 
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where A is the single variational parameter, Q is a func¬ 
tion driven by the many-body potential, and $pc is the 
free Fermi system wave function given by a product of 
Slater determinants, one for each color-flavor combina¬ 
tion of quarks. For a four-body system, the Q function 
takes the following form, depending on the optimal po¬ 
tential at the given configuration, denoted by the sub¬ 
index: 


or 


Q4Q = + h 


3/2 , 3/2 

— -|- r24 , 

(21) 

3/2 , 3/2 

= D4 +^23 

(22) 

3/2 , 3/2 , 3/2 , 3/2 

'2 k +^31 +^4i > 

(23) 


where the 3/2 power is the resemblance of the exact so¬ 
lution of single pairs linked by a linear potential. 

The Slater determinants are composed of single wave 
functions of a particle in a cubic box of side h with Eigen- 
energy 


where nx,ny,nz = 1,2,... We expect this form of the 
variational wave function to be able to describe the de¬ 
parture from the two-mesons case into the four-body 
case, while including the identical particles correlations. 

We can define a particle density parameter as a mea¬ 
sure of the interparticle separation by: p = N/V = 4/L^, 
where N is the number of particles of the system and V 
is the box volume. For a fixed number of particles, the 
change in the particle density corresponds to modify the 
box size and correspondingly to the inter-particle sepa¬ 
ration. 


The first term leads to the noninteracting particle kinetic 
contribution {Tpo) and the remaining terms, after inte¬ 
grating by parts, to: 

( 28 ) 

Thus, using the exponential structure of FA, this term 
corresponds to the expectation value of the following 
function: 




where the sum runs over the coordinates of the four par¬ 
ticles. This term can be seen as the correction to the 
kinetic energy due to the interaction. T hus, the expec¬ 
tation value of the Hamiltonian becomes [23 : 

{H)^ = Tpg + {W)^ + (H);,. (30) 


This structure is useful to find the variational parameter 
A that minimizes it: 


dX 


= 0 . 


(31) 


The free kinetic part is given by the sum of the energies 
of the particles in the box as given by Eqn. ([M]) and is 
independent of A. In order to evaluate <W >, we need 
to know the particles coordinates dependence of Q. For 
the meson configuration, Q is explicitly defined in terms 
of them, while for the tetraquark it involves the vectors 
k and I, whose dependence on the coordinates are not 
analitically known, in this case numerical determination 
is used (see Appendix). 


V. RESULTS 


IV. ENERGY DETERMINATION 


The Hamiltonian of 4 particles of momentum Pi and 
mass rrii interacting through the potential V is: 




(25) 


The expectation values of W and V involves the in¬ 
tegration over 12 variables and can be computed using 
Monte Carlo techniques which, by relying in the impor¬ 
tance sampling, turns the integrals into an evaluation of 
the average of the values of the observable, in particular 
we use a Metropolis algorithm for the sampling, driven 
by the square of the wave function. We computed the 
properties of three systems. The results we present are 
organized as follows: 


Using the above variational wave function, the expecta¬ 
tion value of each term in the kinetic energy can be set 
as: 


• Two mesons. In this case, only meson configura¬ 
tions are allowed in the potential. We characterize 
the meson properties, comparing the results with 
the expected from the analytical solution in the iso¬ 
lated limit and the departure from them as density 
increases. 


where Z(X) = ('I'a|4'a). Defining 'I'a = Fa$, with F\ = 
e~^^, then 

= Fxd^^pG + d^Fx^pG + 2diFxdi^pG. (27) 


• Tetraquark . In this case, only tetraquark configu¬ 
rations are allowed in the potential. We character¬ 
ize the tetraquark properties in the isolated limit 
and the departure from them as density increases. 






5 



FIG. 1: Variational parameter of the two mesons, tetraquark 
and mixed systems as a function of the density, for a. ud = 
uudd mass assignment. They are normalized to the exact 
solution in the isolated meson case Aq. 

• Mixed. In this case, all the possible configurations 
are allowed in the potential. We characterize the 
meson and tetraquark properties and their modifi¬ 
cations due to the presence of the other configura¬ 
tion. 

For all the above cases we explored the mass effect by 
considering the two quarks to be light (u)and the anti¬ 
quarks having variable mass to resemble the md,s,c,b/'m'u 
ratios. Results are presented for several light-quarks to 
heavy-antiquarks mass ratios whenever they are found to 
be relevant. 


A. Variational parameter and energy 

In order to know the best approach of the variational 
wave function to the exact solution, we find the optimal 
value of the variational parameter that minimizes the to¬ 
tal energy at a given density. We repeat the procedure for 
a set of values of the density to know the corresponding 
evolution of the variational parameter. 

In Figure [TJ we show the optimal values of the vari¬ 
ational parameter of the two mesons, tetraquark and 
mixed systems as a function of the density, for a ud = 
uudd mass assignment. The behavior for a light-heavy 
system (ub) is very similar (not shown) thus, for the sake 
of clarity, we only exhibit this case. The results are nor¬ 
malized to the exact solution in the isolated meson case 
Aq. We observe that they exhibit a smooth evolution 
from their asymtotic value at low density decreasing to 
nearly zero at high energy densities, signing the dimin¬ 
ishing effect of the interaction. 

In Figure [2J we have plotted the energy of the sys¬ 
tem {E 4 ), normalized to the exact solution in the iso- 



FIG. 2: Energy of the two mesons, tetraquark and mixed 
systems as a function of the density, normalized to the exact 
solution of two isolated mesons, 2Eo. The upper (lower) panel 
corresponds to a ud (ub) system. 

lated meson case (2i?o) for the three cases: two mesons, 
tetraquark and mixed systems, as a function of the par¬ 
ticle density. The upper (lower) panel corresponds to a 
ud (ub) system. We observe that at low density we re¬ 
produce the exact result for the two mesons system and 
the tetraquark system energy is higher. As the density 
increases, this breach closes and eventually is inverted, 
signing that the tetraquark state is more likely at high 
density. The mixed system turns out to be the one with 
the lowest energy as it takes advantage of both possible 
configurations. 


B. Radial distribution 

Once the optimal variational wave function and energy 
were charaterized, we proceed to compute a set of ob¬ 
servables. In particular, the mean square radius (MSR) 
evolution was obtained, using the optimal pairing infor¬ 
mation. The mean square radius for the meson or the 
tetraquark can be obtained following the general defini¬ 
tion, Eqn. (HU, where the sum runs over the correspond¬ 
ing number of particles. In Figures [3] and 01 we show the 
density evolution of the MSR for mesons and tetraquarks 
respectively, normalized to the corresponding value in the 
isolated case. 

The meson MSR evolution is found to be mass ratio 
dependent. In fact, the light-light combination (ud) de¬ 
creases faster than the light-heavy case (ub), while the 
dependence on being in presence of possible tetraquark 
configurations is mild (solid vs. dashed line) in both 
cases. 

The tetraquark MSR evolution exhibits a mass ratio 
dependence. The light-light system (solid line) increases 
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FIG. 3: Meson mean square radius as a function of density, 
normalized to the isolated case, for light-light (ud) and light- 
heavy (ub) mass ratios. The solid and dashed line correspond 
to the meson and mixed system respectively. 


m2/mi 

2Eo 

EiQo 

(dQ}n 

1 

4.6944 

5.02 ± 0.02 

13.9 ± 0.1 

1.44643 

4.4394 

4.77 ± 0.02 

13.2 ± 0.2 

4.6131 

3.9778 

4.24 ± 0.02 

12.5 ± 0.3 

14.0774 

3.8122 

4.07 ± 0.02 

13.1 ± 0.2 


TABLE II: Two mesons energy, Tetraquark energy and mean 
square radius for different mass ratios respect to the lightest 
u mass, nil in the zero density limit. 


in the low density regime and then rapidly decreases. 
The light-heavy system (dashed line) remains almost un¬ 
changed for a relatively large density region and slowly 
decreases. For both cases we show the modification when 
considering the presence of the mesons (mixed). Both 
become large at low density and rapidly decrease as the 
density increases, approaching the pure tetraquark re¬ 
sult. It is very interesting to find out that during the 
increase of the density, the radial distribution suffers not 
only contractions but also expansions, which are a result 
of a competition between the reduction of the interparti¬ 
cle separation and the decrease of the confining potential 
effect, as it is multiplied by the density decreasing varia¬ 
tional parameter. 

In Table im we show the tetraquark energy {E 4 q^) and 
MSR ((^’Iq)^) in the limit of zero density, for the different 
mass ratios. Twice the meson energy {2Eq) is also exhibit 
as reference. 



FIG. 4: Tetraquark mean square radius as a function of den¬ 
sity, normalized to the isolated case. The solid (dashed) line 
corresponds to a ud (ub) mass ratios, for the Tetraquark and 
mixed systems. 

C. Correlation functions 

Another very useful observable to characterize the 
properties of the system is the two-particle correlation 
function [2^. It measures the probability of finding two 
particles at a relative distance r from each other: 

^ ( E ^ ( 32 ) 

\i<j=i / 

In Figures [5] and |6] we show the meson-meson correla¬ 
tion function at several densities for the two-mesons and 
mixed case, respectively. We exhibit the results corre¬ 
sponding to a light-light (ud) system but a similar be¬ 
havior is observed in the other mass ratio combinations. 
The comparison between both Figures shows that there 
is a modification in the meson-meson correlation function 
by the presence of the tetraquark state at intermediate 
densities, a bump develops in the near tail of the cor¬ 
relation function, driven by the diquark formation (see 
below), and fades out as the density increases and all the 
quark correlations vanish. This may be an alternative 
observable to have indirect evidence of the tetraquark 
formation. Note that at low density the separation be¬ 
tween the two mesons approaches the result obtained by 
Lenz et al 0 of 4.66 rjrM- 

For the two meson case, the quark-quark correlation func¬ 
tion at low density is similar to the meson-meson correla¬ 
tion function (not plotted) as the relevant separation dis¬ 
tance is driven by the separation between mesons, each 
one containing a quark. 















7 



FIG. 5: Meson-meson correlation function at several densities 
(inner labels). Corresponding to a light-light (ltd) system, for 
the two mesons case. 
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FIG. 7: Quark-quark correlation function at several densities 
(inner labels), for a light-light (ltd) system in the tetraquark 
case. 
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FIG. 6: Meson-meson correlation function at several densities 
(inner labels). Corresponding to a light-light (nd) system, for 
the mixed case. 


FIG. 8: Correlation between a quark and all the particles 
regardless of the quantum numbers. Two mesons system. 


In Figure[71 we show the quark-quark correlation func¬ 
tion at several densities (inner labels) for a light-light 
(ltd) system, when only tetraquark configurations are al¬ 
lowed. The diquark formation at low densities is exhib¬ 
ited and, as the density increases, the correlation drops 
signing that, although the tetraquark configuration can 
be identified, the quarks are no longer correlated. Note 
that the diquark size is similar to the bump observed 
in the meson-meson correlation function (Figure 15]) when 
tetraquarks are also formed. At low densities the fraction 
of tetraquarks formed respect to the two mesons is small 


and has no major effect, it is only in the intermediate 
density region where this feature becomes relevant. 

In Figures [5] and [51 we show the correlation between 
a quark and all the particles regardless of the quantum 
numbers, for two mesons and tetraquark systems respec¬ 
tively. Here, we exhibit the case of the light-heavy (uh') 
combination. Although the general features are similar 
for other combinations the specific parameters are not. 
This will be exemplified below, when discussing the zero 
density limit. 

In the zero density limit, the normalized tetraquark 
































FIG. 9: Correlation between a quark and all the particles 
regardless of the quantum numbers. Tetraquark system. 
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FIG. 10: Diquark static structure factor S{q) (a) and correla¬ 
tion function (b) in the zero density limit, for ud (solid line) 
and ub (dashed line) tetraquark systems. 


System 

^0 

Ai 

A 2 

Quark-quark 




ud 

0.64 

1.24 

1.51 

ub 

1.13 

1.1 

1.47 

Quark-all 




ud 

0.24 

1.00 

4.60 

ub 

0.31 

0.98 

3.19 


TABLE III: Diquark and quark-all correlation functions pa¬ 
rameters for ud and wfe-like mass ratios. 


In our case we consider the correlation to depend on the 
magnitude of r. In Figure (TUI we display the diquark 
static structure factor S{q) (a) and correlation function 
(b) in the limit of zero density for both ud (solid line) 
and ub (dashed line) tetraquark systems. 


D. Dynamical recombination 


system correlation function for the diquark (quar-quark) 
and tetraquark (quark-all) can be parameterized by: 

5(’’)g-<? = (33) 

and 

g{r)q-aii = ^o(l + (34) 

In Table Hm we show the value of the parameters of the 
diquark and quark-all correlation functions for the light- 
light and light-heavy cases. Note that the dependence 
on r of the diquark correlation function is close to others 
suggested in the literature [H, |3^ but in our case it is 
dynamically obtained. The short distance region shows 
the competition between the Pauli blocking among the 
identical particles and the attractive potential. 

The static structure factor S'(q) can be obtained as the 
Fourier transform of the correlation function g{r). That 
is 

5(q) = l + yjd^r gir) \ (35) 


The four quark recombination can have important ef¬ 
fects in systems where they are produced very close in 
space. An example can be found in the WW —>■ qqQQ de¬ 
cay, whose spatial separation at LEP2 energies is around 
0.1 fm [T^. Although in the perturbative regime the 
recombination is small, in the non-perturbative regime 
the effect may be important. The typical scenarios to 
estimate the recombination are: Considering spherical 
or elongated bags color sources and the reconnection is 
proportional to the overlap of two color sources; strings 
considered as vortex lines where reconnection takes place 
when the core of the two string pieces cross each other 
M- The recombination between the two mesons con¬ 
figurations is usually implemented [TsI . IT^ . In our case, 
we would have, in addition, recombinations similar to 
tetraquark states which eventually freeze out to two 
mesons. Here, we characterize the modification to the 
two mesons recombination induced by the presence of 
the tetraquark state. 

The model requires that the formed clusters must be col¬ 
orless. The flipping from one configuration to another, 
driven by the minimal potential energy restriction, is a 
measure of the dynamical recombination property of the 
system. We define the recombination probability for the 



























Recombination Probability 
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FIG. 11: Recombination probability as a function of the 
density. The probability for only two mesons recombina¬ 
tion {Pr 2 m.) is exhibited for the two mesons and mixed cases 
(see inner labels). The probability when tetraquarks are ac¬ 
counted as a recombination stage (Prmix) corresponds to the 
mixed case. The solid and dashed lines are for ud and ub 
systems respectively. 


two mesons case by: 

N{Vml O V^2) 
N{Vml)PN{Vm2) 


(36) 


where N{Vi —> Vj) denotes the number of flippings from 
the i to the j configuration, and N{Vi) is the number of 
times the system visits the i configuration. If we consider 
the tetraquark configuration as another possible configu¬ 
ration contributing to the recombination, the probability 
becomes: 


Pfmix — {N{Vml ^ Vra2) + N{Vml ^ V4 q) + 
N{Vm2 ^ ^ q ) + N{V4Q —>■ Vaq)] / 
[N{Vml) + N{Vm2) + N{V4 q)] (37) 

In Figure HH we show the probabilities as a function 
of the density. The probability for only two mesons 
recombination {Pr 2 m) is exhibited for the two mesons 
and mixed cases (see inner labels). The probability 
when tetraquarks are accounted as a recombination stage 
{PTmix) corresponds to the mixed case. The solid and 
dashed lines are for ud and ub systems respectively. 

A qualitative estimate of this effect can be exhibited 
considering a simple expansion model, such that we can 
evolve the system along the density profile. We can define 
the color strength function by: 

n(a;) = Pfrag{t)Pr{x) (38) 

where Pfrag{t) = ^xpi—t^/Tfrag) i® probability that 
the system has not yet fragmented, with Tfrag ~ 3r?i the 


FIG. 12: Golor strength due to the recombination probability 
in the two mesons and mixed cases as a function of the radial 
size of the system (in units of the meson radius). The solid 
and dashed lines are for ud and ub systems respectively. 


proper lifetime, taken as three times the meson radius 
rfi [i^. The recombination probability Pr{p) can be set 
in terms of the radius of the sphere (a;) corresponding to 
the given density p = . In Figure [T^ we plot 

the evolution of the color strength as a function of the 
radial size of the system, in units of the meson radius. 
The case when considering the tetraquark and meson re¬ 
combination [Pumix) largely modifies the estimate when 
the tetraquark state is not included {Pr 2 m)- The flavor 
composition modifies the behavior at intermediate dis¬ 
tances, and the total effect becomes negligible after four 
times the hadron radius, corresponding to nearly the size 
of the end of the overlap of two hadrons. 


E. Tetraquark potential 


The tetraquark potential depends on the quarks posi¬ 
tions and two auxiliary vectors, placed in such away that 
the total length linking the quarks is the shortest one. 
These auxiliary vectors are modified in a non-trivial way 
whenever a single quark changes its position. An effec¬ 
tive behavior of the potential can be set as linear respect 
to the invariant length, R = ’ 


,2 . 


V{R) = Ro + BR 


(39) 


where i?o is the value of the potential at zero distance, 
which is expected to be modified by the short distance 
coulomb-like correction In Figure [131 show the 
slope B of the linear behavior from the simulation as a 
function of density for the case when only tetraquarks 
are allowed to form (solid line) and when mesons and 
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FIG. 13: Slope (B) of the effective tetraquark linear potential 
as a function of density, for the pure tetraquar system (solid 
line) and meson- tetraquark mixed system (dashed line). 

tetraquarks are allowed to form (dashed line). The slope, 
in the tetraquark case shows small dependence on the 
density, and B{p « 0) = 0.84 ± 0.02. In the mixed case 
there is a significant density dependence, starting below 
the corresponding value at zero density {B{p Ri 0) = 
0.67 ± 0.02) and as the density increases they approach 
to each other. This behavior is similar when considering 
different flavor systems. 

In Figure [HI we show the effective four-body contact 
potential, V 4 (contact), as a function of density, deter¬ 
mined as the average potential irrespective of the length 
linking the particles. The pure tetraquar system (solid 
and dashed lines for ub and ud flavors respectively) ex¬ 
hibit a dependence on the mass ratio in the low energy 
regime and, as density increases, they become similar. 
The meson-tetraquark mixed system (dotted line) has no 
dependence on the mass ratio (labeled as mixed), but is 
is very sensitive to the presence of the mesons, becoming 
very large in the low energy regime. This behavior cor¬ 
responds to the increasing distance among quarks which 
are influenced by the two mesons configuration that is 
trying to bring them apart. 


VI. CONCLUSIONS 

We have performed a Monte Carlo simulation of two 
identical light quarks qq and two identical antiquarks 
QQ, considering three possible structures: two mesons, 
tetraquark and mixed configurations. We determined 
whether it is energetically more favorable to form a 
tetraquark or two mesons and the mixing among them, 


FIG. 14: Effective tetraquark contact potential as a function 
of density, for the pure tetraquar system (solid and dashed 
lines for ub and ud flavors respectively) and meson- tetraquark 
mixed system (dotted line). 


as a function of the particle density. Several light-quarks 
to antiquarks mass ratios were used and exhibited when¬ 
ever they were found to be relevant. The meson MSR 
evolution was found to be mass ratio dependent. In fact, 
the light-light combination {ud) decreases faster than the 
light-heavy case {ub), while the dependence on being in 
presence of possible tetraquark configurations is mild. 
The tetraquark MSR evolution also exhibits a mass ra¬ 
tio dependence, and turns out to be very sensitive to the 
presence of the mesons. 

We have shown that there is a modification in the 
meson-meson correlation function by the presence of the 
tetraquark state at intermediate densities, where a bump 
develops in the near tail of the correlation function, 
driven by the diquark formation. To gain insight into 
the features of the tetraquark state, we computed the 
quark-quark correlation function at several densities. We 
observed the diquark formation at low densities and the 
drop on the correlation as the density increased. The 
quantum correlations among identical particles appeared 
as Pauli blocking in the short distance regime. In the zero 
density limit, we determined the correlation functions be¬ 
tween quarks and the correlation between a quark and 
the rest of the particles. A parameterization was found, 
which is useful to compute additional static properties, 
in particular we computed the diquark static structure 
factor, which exhibited a small mass ratio dependence. 

We did track the dynamical flipping among configura¬ 
tions and determined the recombination probability evo¬ 
lution as a function of the particle density. We have 
shown that the probability is largely affected when con- 
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sidering the tetraquark as an intermediate recombination 
state. Using a simple model evolution of the system along 
the density profile we have shown that the recombination 
effects remain important beyond a half the overlap of two 
mesons. 

The linear behavior of the four-body potential on the 
invariant length linking the quarks was analyzed and 
found that the presence of a mixed state is reflected in the 
decreasing strength of the slope. We determined an effec¬ 
tive four-body contact potential as the average strength, 
irrespective of the length linking the particles. This is 
very sensitive to the presence of the mesons and when 
they are not considered there is a dependence on the 
mass ratio in the low energy regime. 

The approach we have followed to the tetraquark sys¬ 
tem is intended to exhibit the gross features. Certainly, 
the short distance information and more elaborated in¬ 
teractions should have an effect on the system energy and 
the corresponding spectrum. In this work, we have lim¬ 
ited ourselves to the case of identical quarks and identical 
antiquarks. Different particles effects should be reflected 
in the high density regime due to the absence of Fermi 
correlations. 


Appendix: Kinetic term 


Analogous results are obtained for the other particles, 
adding all them up we get: 


{W{Qml)) \ — 


9A^ 


f <ri3> 
V M13 


< ^24 > \ 
M24 ) 


(A.l) 


where is the corresponding reduced mass. The Qm^ 
configuration kinetic term is obtained by the 3 0 4 ex¬ 
change. 

Tetraquark: In this case we have two additional vectors 
k and I, whose explicit dependence on the coordinates of 
the quarks and antiquarks are unknown. Let us illustrate 
this implication in the contribution from the coordinate 
X of particle I: 


(rjf -t -f r^^ + -b 


— n t -b 

~ l' Ik ’-’xi'lk -T ... 


2' Ik 

The first term derivative can be set as: 


(A.2) 


^xi^lk — [(^1 ^k^ (1 ^x\^k) iyi Uk^ ^xiVk 

-{zi - Zk)dxiZk\ /rIk (A.3) 


The kinetic term W Eqn. (uni) is evaluated as follows: 

Meson: In this case we have two configurations, corre¬ 
sponding to Qmi and Qm 2 ■ Let us illustrate the form of 
the contribution from the coordinate x of particle I, by 
taking the derivative of the first element of Qmi '■ 


where we do not know the derivative of the k coordinates 
respect to the xi coordinate. The derivatives are numer¬ 
ically evaluated as dx^Xk ^ Axfc/Axi. Analogous terms 
are evaluated for all the rest of the particles coordinates. 


,, „3/2 _ 3 (xi - Xg) 

0'a:i7’i3 - 2 1/2 ■ 

D3 

Thus, if we compute the contribution for the 3 coordi¬ 
nates and add their squares, we get: 


[dxiQrmf + [dy^Qrmf + [dziQrmf = 
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